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The nonlinear evolution of resistive double tearing modes (DTMs) with safety factor values q = 1 
and g = 3 is studied with a reduced cylindrical model of a tokamak plasma. We focus on cases where 
the resonant surfaces are a small distance apart. Recent numerical studies have shown that in such 
configurations high-m modes are strongly unstable and may peak around m = rripeak ~ 10. In this 
paper, it is first demonstrated that this result agrees with existing linear theory for DTMs. Based 
on this theory, a semi-empirical formula for the dependence of mpeak on the system parameters is 
proposed. Second, with the use of nonlinear simulations, it is shown that the presence of fast growing 
high-m modes leads to a rapid turbulent collapse in an annular region, where small magnetic island 
structures form. Furthermore, consideration is given to the evolution of low-m modes, in particular 
the global m = 1 internal kink, which can undergo nonlinear driving through coupling to fast 
growing linear high-m DTMs. Factors influencing the details of the dynamics are discussed. These 
results may be relevant to the understanding of the magnetohydrodynamic (MHD) activity near 
the minimum of q and may thus be of interest for studies on stability and confinement of advanced 
tokamaks. 



I. INTRODUCTION 

The advanced tokamak (AT) scenario, where the max- 
imum current density is located off the magnetic axis, is 
of considerable interest for achieving thermonuclear fu- 
sion conditions and a quasi-steady-state operation in fu- 
ture tokamak devices (e.g., Refs. P, The associated 
non-monotonic q profiles have pairs of resonant surfaces 
with the same rational value qs = mjn and may give rise 
to double tearing modes (DTMs) Q. Several detailed 
studies of DTMs were motivated by their possible role 
in rapid current penetration, compound sawtooth oscil- 
lations, off-axis sawtooth crashes and disruptions (e.g., 
Ref. and references therein). It is important to un- 
derstand the behavior of these instabilities in order to 
ensure efficient profile control and safe operation of toka- 
mak devices. 

In this paper we consider cases with q^ = 1 and qs > 1 
DTMs, both of which have attracted much attention in 
experiments. For instance, compound sawtooth oscilla- 
tions were observed in Tokamak Experiment for Tech- 
nology Oriented Research (TEXTOR) after two gs = 1 
resonances had formed Off-axis sawteeth were ob- 
served on Tokamak Fusion Test Reactor (TFTR) when 
the minimum of the safety factor, gmin, is near or below 2 
. Resistive and neoclassical gs = 2 DTMs were investi- 
gated in Axisymmetric Divertor Experiment (ASDEX)- 
Upgrade, one motivation being the possible interaction 
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of these modes with internal transport barriers (ITBs) 
0, In Japan Atomic Energy Agency Tokamak Up- 
grade (JT-60U) qs 3 DTMs are thought to play a 
crucial role in disruptions, as indicated b y ex perimental 
observations and numerical results [§, [lol. Illl| . 

In contrast to previous numerical studies, which fo- 
cused mostly on cases where two resonant surfaces are lo- 
cated a relatively large distance apart (e.g., Ref. [ll]) and 
the linearly most unstable mode has the lowest possible 
poloidal mode number [e.g., (m, n) = (3,1) for = 3], 
our interest lies in the regime where the distance between 
neighboring resonant surfaces is still small. For this case, 
it was found that DTMs (and multiple tearing modes in 
general) with high poloidal and toroidal mode numbers m 
and n are strongly unstable . The linear instability of 
DTMs for equilibria with small inter-resonance distances 
was studied numerically in Ref. [J- The purpose of the 
present paper is (i) to shown that the findings of Ref. [3| 
agree with linear DTM theory, and (ii) to present first 
nonlinear simulation results involving high-m DTMs for 
cases with qs = 1 and > 1- For simplicity, a reduced 
magnetohydrodynamic (RMHD) model is employed. 

Existing linear theory for DTMs [l3| predicts that in 
the strongly coupled limit (small mode numbers, small 
inter-resonance distance) the linear growth rate increases 
with the mode number as 7iin c>c m^l^ . Here, it is 
shown that this scaling agrees with numerical results. 
A semi-empirical analytical formula for estimating the 
mode number mpeak of the fastest growing mode in terms 
of the system parameters is proposed, based on the tran- 
sition criterion between the strongly and weakly coupled 
limits derived in Ref. [l3| . 

The magnetic island dynamics during the nonlinear 
evolution of (js = 3 and q^ = \ DTMs is described. In 
the (js = 1 case, the m ~ I internal kink mode is unsta- 
ble and eventually dominates the dynamics. However, 
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full reconnection is not possible for small inter-resonance 
distances, so the core merely undergoes oscillatory mo- 
tion. Full reconnection is only possible when the inter- 
resonance distance has already increased beyond the limit 
where modes with to > 1 become sub-dominant or negli- 
gible compared to the m = 1 mode. This is demonstrated 
using an intermediate case where the m = 2 mode has a 
slightly higher growth rate than the to = 1 mode. 

Finally, the early evolution of the to 1 mode is inves- 
tigated in detail. Of particular interest is the nonlinear 
driving due to fast growing high-m DTMs. It is shown 
how the efficiency of this driving depends on the initial 
conditions of the simulation and that the driving, despite 
its radial localization, is capable of triggering the global 
resistive m = 1 internal kink mode, similarly to the case 
of Q s = 1 triple tearing modes (TTMs) studied recently 

This paper is organized as follows. In Section |TT] the 
physical model is introduced. In Section Hill we describe 
the equilibrium configurations used and their linear dis- 
persion relation. Section IIVI is dedicated to a compar- 
ison between linear theory and numerical data, and in 
Section |V] nonlinear simulation results are presented. In 
Section rvTl we draw conclusions, discuss possible applica- 
tions and motivate further research in this direction. 



II. MODEL 

We use the reduced magnetohydrodynamic (RMHD) 
equations in cylindrical geometry and in the limit of zero 
beta [H, • This model has proven to be useful in stud- 
ies of MHD instabilities, when the focus is on a qualita- 
tive description of fundamental aspects of the magnetized 
plasma system, as is the case here. The RMHD model 
governs the evolution of the magnetic flux function -0 and 
the electrostatic potential 0, as described previously in 
Ref. [J. The normalized RMHD equations are 

= [^,^]-dc^-S^lifij-E„) (1) 
dtu = [w, 0] + [j, iP] + da + Re^lVlu. (2) 

The time is measured in units of the poloidal Alfven time 
''Hp = \/ l^-oPmO-/ Bq and the radial coordinate is normal- 
ized by the minor radius a of the plasma. p,„ is the mass 
density and Bq the strong axial magnetic field. The cur- 
rent density j and the vorticity u are related to '0 and (f> 
through j = — V^V and u = 'S/\4', respectively. 

Resistive diffusion is measured by the magnetic 
Reynolds number S'jjp = T,,/rHp in Eq. ([T|). with t,, = 
a^/zo/?/o being the resistive diffusion time and 770 = 'rj(r ~ 
0) the electrical resistivity at r = 0. We use S'np = 10^, 
which is numerically efficient and physically reasonable 
in the framework of the model used. Viscous dissi- 
pation is measured by the kinematic Reynolds number 
R^Hp = ct^/i'THp in Eq. where v is the kinematic ion 
viscosity. We choose regimes where the Prandtl number 



Pr = S'Hp/-ReHp satisfies Pr ~ 10~^, so that the vis- 
cosity effect is limited to small-scale flows and does not 
affect the instability of the dominant modes. 

The source term S^^Ea in Eq. ^ compensates the re- 
sistive diffusion of the equilibrium current profile on the 
time scale tr/thp = -D^2'^Hp, where D12 is the length 
scale of interest, i.e., here the inter-resonance distance 
(normalized by a). Eq is taken to be constant, so the re- 
sistivity profile is given in terms of the equilibrium cur- 
rent density distribution as 77 (r) = j{r ~ 0)/j(r). For 
simplicity, the temporal variation of the resistivity pro- 
file fj is neglected. 

As in Ref. Q], each field variable / is decomposed into 
an equilibrium part / and a perturbation / as 



fir,d,C,t)^ f{r) + fir,^X,t). 



(3) 



The system is described in terms of the Fourier modes, 
iprn.n and (f>m.m obtained from the expansion 

fir, §X,t)^^J2 + c.c, (4) 

with TO being the poloidal mode number and n the 
toroidal mode number. In the following, the (to, n) sub- 
scripts will often be omitted for convenience. We consider 
only the nonlinear couplings between modes of a single 
helicity h = m/n, so the problem is reduced to two di- 
mensions. Results for the linearized system are obtained 
usingj_nitial-value and eigenvalue solvers as described in 
Ref. The nonlinear RMHD equations are solved nu- 
merically using the simulation code described in Ref. . 

In order to ensure numerical accuracy, we have evalu- 
ated the energy balance (temporal change in the system's 
energy compared to the dissipated energy) and compared 
the results obtained with different numbers of grid points 
and Fourier modes. In particular, the time histories and 
mode structures of individual modes were inspected in 
detail. Both linear and nonlinear calculations were per- 
formed with a grid spacing of Ar = 5 x 10~* and smaller. 
The number of Fourier modes is specified below for each 
case. In the regimes where the linear theory of DTMs 
[1^ is valid, it was used to benchmark numerical re- 
sults (e.g., the resistivity scahng S'gp with a e [1/3, 3/5]; 
cf. Figs. 11-14 in Ref. %. 



III. EQUILIBRIUM AND LINEAR 
INSTABILITY 



The equilibrium state is taken to be axisymmetric 
(only m = n = components) and free of flows, i.e., 
(j) ~ u ^ 0. The equilibrium magnetic configuration is 
uniquely defined in terms of the safety factor q{r), and 
the magnetic flux function and current density profiles 
are obtained though the relations 



q 



= — ^-V'o.o and jo,o = -^j ■ (5) 

r ar r ar q 
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The model equation used for the q profile is 



(7(r)=qo.J^i(r).{l + (r/ro)''"('')} 



l/w(r) 



(6) 



where 



w{r) 
Fi{r) 



= f A 



-l/[2w(rA)] 



Wo + Wir , 

l + /iexp{-[(r-rn)/ri2]'} 



With the parameter values in Table U the equilibrium q 
profiles shown Fig. [T] are produced, each of which has 
two resonant surfaces with qs = q{rsi) = m/n at the 
radii r ~ r^i {i = The distance between the reso- 

nances, Di2 ~ |r's2 — Tgil was chosen sufficiently small, so 
that broad spectra of DTMs are unstable, with dominant 
modes having m greater than the lowest poloidal mode 
number that is consistent with a given field line pitch 
qs = m/n. The dispersion relations (spectra of linear 
growth rates) ^iinim) are plotted in Fig. [21 In Cases (D- 
1) and (D-2), pairs of resonant surfaces with qs — i and 
1, respectively, are located a small distance D12 — 0.06 
apart, so that the fastest growing mode has mpoak = 8. 
In addition, we consider Case (D-3) with qRes = 1, which 
has a larger inter-resonance distance, D12 = 0.21. Here, 
the dominant mode is mpoak = 2. The characteristics 
of all cases are summarized in Table HIl Linear eigen- 
mode structures of DTMs with various mode numbers m 
were presented in Ref. [1]. Note that there can be up to 
two unstable eigenmodes for a given (to, n). Eigcnmodes 
peaking in the region < r < Tgi (for qs = 1) or near rgi 
(for (js > 1) are denoted by M'^^^ while those extending 
to the outer resonant surface are labeled M^'^\ When 
D12 is small, the Af'^^^-type eigcnmodes are usually dom- 
inant for TO > 2 and only their growth rates are plotted 
in Fig. H 

Finally, a comment is in place with regard to the re- 
gion that can undergo magnetic reconnection in the cases 
studied here. The profiles of the equilibrium helical fiux 
functions V'* = + r^/i'^qs) m Fig. [3] show that the re- 
connectable regions do not include the magnetic axis at 
r = 0. Since all cases satisfy V'*(z_= 0) < V'*(^s2), the 
reconnection can only be partial [isl . [Toj . 



IV. COMPARISON WITH LINEAR THEORY 

Pritchett. Lee and Drake (PLD) developed a lin- 
ear theory for resistive DTMs by applying the techniques 
and results from the theory for the resistive m ~ 1 in- 
ternal kink mode by Coppi et al. [13]. PLD derived an- 
alytical forms for the linear growth rate 7iin of DTMs in 
the limit of strong and weak coupling. The strongly cou- 
pled limit applies to modes which are close to idcal-MHD 
marginal stability, whereas weakly coupled modes are 
strongly idcal-MHD stable. Defining the poloidal wave 



Case (D-2) 

/ Case (D-3) 




FIG. 1: (Color online). Equilibrium safety factor profiles q{r) 
for Cases (D-l)-(D-3). The properties of these profiles are 
listed in Table lU and dispersion relations are shown in Fig. [21 




poloidal mode number m 

FIG. 2: (Color online). Spectra 7iin(?ri) of unstable DTMs 
for the three cases in Fig. [H for Shp = 10® and Reup — 10* 
[i?eHp = 10'' in Case (D-3)]. For Case (D-2) (circles) the 
growth rates of both m — 1 eigenmodes, M'^' (single kink) 
and M'^' (double kink), are shown: ^'•^''{m = 1) = 4.0 x 10"^ 
and 'y^^\m = 1) = 2.5 x 10~^. All other growth rates belong 
to M^^-'-type eigenmodes, as defined in the text. 



number fc^ = m/ro [tq = (rgi + rs2)/2] and the distance 
Xs = D12/2, the validity of the dispersion relations in the 
two limits is given by the constraint [l^ 




^ k^Xs ^ 1 



weak 



(7) 



where B' 



s(rsi)/qs, with s — rq'/q being the mag- 
netic shear evaluated at a resonant surface. The lower 
limit is based on the requirement that the resistive layer 
width Sri be smaller than the distance Xs, and the upper 
limit corresponds to a largc-aspect-ratio approximation 
used in the derivation of the dispersion relation. Modes 
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TABLE I: Parameter values for the q profiles in Fig. [T] using the model formula ((6}. 



Case 


Qo 




Wo 


Wi 


m 


n 


fi 


rii 


ri2 


(D-1) 


2.6 


0.655 


3.8824 





3 


1 


-0.238 


0.4286 


0.304 


(D-2) 


1.3 


0.655 


3.8824 





1 


1 


-0.238 


0.4286 


0.304 


(D-3) 


1.25 


0.655 


3.8824 





1 


1 


-0.238 


0.4286 


0.304 



TABLE IL Properties of the q profiles shown in Fig. [T] The linear instability characteristics of these cases were previously 
studied in Ref [J| , with the case labels given in the second column. The values of the magnetic shear and the resistivity at the 
resonant surfaces are denoted by Si = s{rsi) and r}i = fj{rsi), respectively. The mode numbers of the fastest growing mode, 
'Tipcak, are valid for Sup ~ 10^ and Reap = 10* [-RcHp = lO'^ in Case (D-3)]. 



Case 


Case in Ref. [4] 




Qmin 


Di2 


Sl 


S2 


Vi 


m 




(D-1) 


(nib) 


3 


1.99 


0.06 


-0.10 


0.12 


0.76 


0.84 


8 


(D-2) 


(la) 


1 


0.99 


0.06 


-0.10 


0.12 


0.76 


0.84 


8 


(D-3) 


(lb) 


1 


0.96 


0.21 


-0.20 


0.45 


0.75 


1.07 


2 



D,, = 0.21 ,' / 




(D-3)' / 




>J / 


~ 




V*(r=0) 


^*(^S2) / 




/(D-1), (D-2 





M 1 






Cases (D-1), (D-2) 


D =0.06 


Case (D-3) 





0.4 0.5 

r 



FIG. 3: (Color online). Equilibrium helical fiux functions 
'0j(r) for Cases (D-2) and (D-3). The radial extent of "recon- 
nectable" regions is indicated by arrows. 



that fall into the transitional regime between the strongly 
and weakly coupled limits are neither close to marginal 
stability nor strongly ideal-MHD stable. In this section 
we reproduce essential steps from PLD's derivation and 
compare the theoretical predictions with numerical re- 
sults, focusing in particular on the m dependence of the 
linear growth rate. It is then shown that a semi-empirical 
formula for the mode number mpoak of the fastest grow- 
ing mode can be extracted from the theory. 

The dispersion relation for strongly coupled DTMs 
(small fc^Xs) is [13] 



7ii. 



(strong couphng), (8) 



where A = A(Ah) is a normalized growth rate. The de- 
pendence of A on the ideal-MHD instability drive mea- 
sured by Ah was obtained by PLD through asymptotic 



matching. At marginal stability one has Ah = and 
A(0) = 1. Equation ([5]) implies that the growth rate in- 
creases with the poloidal mode number m. A comparison 
of the rr?l^ power law with numerically computed growth 
rate spectra from Ref. [3] is given in Fig. [H Good agree- 
ment is found, apart from the data point 7hn(TO = 1) 
in Case (D-2) with ^Hp = 10^ [arrows in Figs. U] (a) and 
(b)]. This deviation is most likely related to the fact that 
(5^ = YixJ{klB;^SYL^)fl^ w 0.37x1)12 = 0.74xx, in that 
case, which violates the condition 5^ ^ [13] • Here, (5^ 
is the linear resistive layer width given by Eq. (22) in 
Ref. [3]- Note that in order to fit Eq. ^ to the numeri- 
cal data in Fig. 2] we use A < 1, which means that DTMs 
are only approximately marginally stabile and lie slightly 
in the ideal-MHD stable domain. 

In the limit of weak coupling (large fc^Xg) the disper- 
sion relation is VA 



7iii 



V7hr(-i/4); 

The quantity 7h is given by 



/ 8r(5/4) ^fc^^f x3/5 



5, 



Hp 



7h 



Jdx.i?.^(x), 



(9) 



(10) 



and thus depends on the shape of the q profile [here, 
= B^{r) — B^{rsi) is the helical field that reverses 
sign across r = rg,]. In the limit of Xg — *■ one can 
approximate by a parabola centered half-way between 
the resonances, w B'^{^^ ~ a;^)/2xs, which yields 7h ~ 
B'^klxl- With this, Eq. ^ gives 



7lin oc m 



-6/5 



(weak coupling), 



(11) 



[Note that the parabolic approximation affects only the 
Xs dependence, not the m dependence. Here, we assume 
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that Xs be small and the weak coupling is realized through 
large to.] This result shows that in the weakly coupled 
limit the growth rate of DTMs decreases with increasing 
TO. The m~^/^ power law does not fit the data, which 
can already be seen from the convexity (d^7iin/dTO^ < 0) 
of the spectra for the cases (D-1) and (D-2) in Fig. [2] as 
opposed to the concavity of 7iin in Eq. (fTT|) . Here, it is 
useful to check the validity of the weakly coupled limit. 
Substitution of the profile parameters Xg ~ 0.03, vq ~ 

0. 42 and « 0.11 into Eq. we find 7 < TO„cak < 14 
(Am < 7) for ^Hp = 10^. The fact that the to^^/s law 
is not observed in our data suggests that the transition 
regimes associated with the lower and upper bounds of 
TOweak are so broad that they overlap for Am < 7 (note 
that the theory docs not define the exact meaning of the 
symbol "^"). For larger values of S'np the lower limit for 
TO is reduced. Indeed, the high-TO side of the spectrum for 
the cases (D-1) and_(D-2) becomes concave for S'np = 10* 
(cf. Fig. 10 in Ref. where 5 ^ mwoak ^ 14 is required 
according to Eq. ([7]). For case (D-3) where Xs « 0.1, 
the weak coupling limit is valid for 1 <C TOwcak ^ 4 
(•^Hp = 10^) which encompasses essentially all unstable 
modes as can be seen in Fig. [5] Again, the number of 
unstable modes is too small for the spectrum to follow the 
TO~^/^ law in any significant interval. We conjecture that 
an appreciable number of modes satisfying the dispersion 
relation Eq. ^ for the weakly coupled limit may only be 
found for very small Xs and very large S'np. 

One more comment is in place with regard to the upper 
boundary of the weak coupling limit, k^Xs 1. When 
this boundary is exceeded by increasing Xg then weak cou- 
pling implies that the modes are spatially decoupled and 
eventually turn into independent single tearing modes. 
On the other hand, if a;s is left sufficiently small, such 
that the perturbation at one resonant surface can still be 
"felt" at the other resonance, an increase in = m/ro 
such that ro/m < Xg may be expected to weaken the 
mutual driving of the perturbations: an 0-point at one 
resonance will "push" not only the facing X-point but 
also on 0-points. To the best of the authors' knowledge, 
this regime is not yet understood. The stabilization of 
high-TO modes seems to depend on other factors, includ- 
ing the mechanism breaking the ideal-MHD constraint. 
For instance, coUisionlcss DTMs due to electron inertia 
tend to have a broader spectrum of unstable modes than 
resistive DTMs for the same q profile plj . 

The growth rates plotted in Fig. [2] increase with to 
for TO < TOpoak and a decrease for to > TOpoak- Ac- 
cording to Eqs. ([H]) and (fTTj) . this behavior corresponds 
to that predicted for the strongly and weakly coupled 
limit, respectively. Despite the lack of a distinguished 
band of weakly coupled modes in the cases considered 
here, we conjecture that an estimate for the mode num- 
ber of the fastest growing mode TOpoak may be obtained 
from the transition criterion in Eq. ^ derived by PLD, 

1. e., k^Xs ~ [k'^/{B'^Siip)Y^^ ■ Solving this relation for to 
yields 

TOtrans « To / (x^ B',Snpy ■ (12) 




FIG. 4: (Color online). Comparison between theory and sim- 
ulation: m dependence of the linear DTM growth rate. Only 
the growth rates of Af'^'-type modes are shown. The dashed 
lines indicate the scaling law 7iin oc m^^^. To fit Eq. ((Sjl to 
the data, A in the range 0.7 < Ah < 0.9 is used. 

(a) : Spectra 7iii'(m,5Hp) of Case (D-2) for S'np = lO'', 10^ 
10*. The data are the same as in Fig. 10 of Ref. 

(b) : Spectra 7iin (™-i'?s) for Qb = 1 and 3 [Cases (D-2) and 
(D-1)] with Shp = 10^. The data are the same as in Fig. 6(a) 
of Ref. [|. 



for mode numbers in the transitional regime. Let us com- 
pare some values obtained from Eq. (|12p with numerical 
results. For instance, for Case (D-2) (where 7'o = 0.42, 
Bg « 0.11, Xs w 0.03) one obtains TOtrans — 7,5,4 for 
Siip = 10^ 10^ 10*. The measured values for 

^peak are 

8, 6, 5 (Fig.m see also Fig. 10 in Ref. [3]), which suggests 
that TOpeak ~ 'Titrans + 1 ■ Tests with Other configurations 
gave similarly good agreement, despite the fact that the 
Xs dependence is described only approximately under the 
assumption that g(r) is parabolic around r ~ rQ. Indeed, 
there is not much freedom for varying the shape of the 
q profile in the inter-resonance region when the distance 
between the resonances must be small, so for realistic q 
profiles the parabolic approximation may be expected to 
be sufficiently accurate. 

Based on the good agreement between linear theory 
and simulation, we propose the semi-empirical formula 

"^-^^ " ix!B^L,V/^ + ^'^^ 

for the dependence of TOpeak on the system parameters. 
Equation ()13p is useful for small D12 = 2xs, where 
"^poak > 1- Note that a small inter- resonance distance 
also implies that the difference between the magnetic 
shears si and S2 is small. Hence, due to the weak shear 
dependence it is not so important whether TOpoak is eval- 
uated using si or S2- 

The poloidal mode number of the fastest growing lin- 
ear mode, TOpoak, is useful for the interpretation of the 
nonlinear dynamics, since it determines the size of the 
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magnetic island structures. Numerically computed val- 
ues for the profiles used in this study are given in Table [III 



V. NONLINEAR RESULTS 

Starting from an unstable equilibrium, the instability 
growth is excited by applying an initial perturbation of 
the form 

= 0) = J E - l)e^i---^'+^'>-) + c.c, (14) 



where is the perturbation amplitude (typically 

10^''), -d^ = -d — Qs^C is a helical angle coordinate and 
•dom is an initial phase shift. The values i?om = and tt 
are assigned to each m in a random manner. This intro- 
duces some degree of incoherence while retaining mirror 
symmetry about the x axis (due to parity conservation in 
RMHD). This restriction improves numerical accuracy, 
simplifies visualization and has no significant effect on 
the central claims of this paper. 

We begin with a description of the magnetic reconnec- 
tion dynamics in Sec. IV Al where the system as a whole is 
considered. In Sec. lVBl the evolution of individual modes 
is analyzed in detail. 



A. Magnetic reconnection dynamics 

First, consider Case (D-1) where two (/s = 3 resonances 
are located a small distance Z?i2 = 0.06 apart. A non- 
linear simulation was carried out using \E'o,m>o = 10~^, 
Sup = 10^ and i?eHp = 10^, and including 32 modes 
(m = 0,3,..., 93). Snapshots showing contour plots of 
the helical flux function -0, = i/' + ?'^/(2gs) and the elec- 
trostatic potential (j), as well as instantaneous profiles of 
q{r) and (j*) = [-V5^V*]o,o = jo.o - 2/gs, are presented 
in Fig. [5] We observe that, in response to the random 
broad-band perturbation applied, magnetic reconnection 
occurs simultaneously at many locations, giving rise to a 
multitude of small magnetic islands. In Figs. EJA) and 
(B) it can be seen that the dominant island sizes cor- 
respond to the mode numbers m = 6 and 9 (dispersion 
relation: mpoak ~ 9). The onset of the reconnection is 
determined by the growth rate of the (9, 3) mode, which 
is about 1.6 times that of the (3, 1) mode (cf. Fig. 
The reconnection leaves behind an annularly flattened q 
profile, as can be seen in Fig.[5]JC). 

Next, let us investigate the response oi Qs — 1 DTMs 
in Case (D-2). The calculation was performed with 
^'o,m>o = 10"^, S'hp = 10^ and i?eHp = 10^, and 128 
modes (m = 0, 1, 127) were included. As can be seen 
in Fig. [BfA), the reconnection dynamics begin with an 
annular collapse with more or less turbulent patterns 
dominated by mode numbers around m ~ 7-9, in ac- 
cordance with the peak of the linear dispersion relation 
(wpoak = 8). Note in Fig. EI A) that the q profile has 
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FIG. 5: (Color online). Reconnection dynamics with gs = 3 
DTMs for small inter-resonance distance D12 = 0.06 [Case 
(D-1)]. The snapshots were taken at (A) t = 1300, (B) t = 
2000 and (C) t — 2500. Each snapshot consists of contour 
plots of the helical flux V'* (top) and the electrostatic potential 
(j> (bottom). Arrows indicate the flow directions. On the 
right-hand side, the instantaneous profiles q{r,t) and (j*) = 
[j, (r, t)]o,o are shown. 5hp ~ 10^, i?eHp = 10*. 
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been flattened in the inter-resonance region and that the 
resistive m = 1 internal kink mode is not yet involved in 
the dynamics. The kink appears at a later time, leading 
to a growing core displacement inside the turbulent re- 
gion, as can be seen in Fig. [HlB). After all reconnectable 
flux surfaces have been reconnected a rebound occurs and 



the core displacement decays as indicated by the arrows 
in Fig.[6l^C). The further evolution was not investigated, 
but continuing oscillation of the core is to be expected. 

Finally, in Fig. [7] simulation results are presented 
for Case (D-3). The calculation was performed with 
*o,m>o = 10"'^, 5'hp = W and i?eHp = lO'^, including 
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128 modes (m = 0, 1, 127). Case (D-3) is a realization 
of the intermediate regime where two = 1 resonances 
are located so far apart that the fastest growing mode is 
7Ti = 2, closely followed by ?ti = 1. It can be seen that 
magnetic islands with dominant mode numbers m = 2 
and m = 1 largely determine the structure of the mag- 
netic surfaces [Fig.[7]^A) and (B)]. In the present case, the 
m = 1 mode tends to dominate on the outer resonance 
and m = 2 on the inner one, which leads to a strong 
deformation of the inter-resonance region into a D-shape 
[Fig. EfC)]. This calculation had to be terminated soon 
after snapshot (C) due to a continuing increase in the 
energies of high-m modes. An island separatrix is ap- 
proaching the coordinate origin, r = 0, and our numer- 
ical code is not suitable for further following dynamics 
with such a kind of asymmetry. 

Equilibria with larger inter-resonancc distance D12 
than in Case (D-3) have a dominant m = 1 mode and 
the dynamics proceed as described in earlier studies (e.g., 
Refs. i,[i3). 



B. Detailed evolution of individual modes 

For the discussion of the evolution of individual Fourier 
modes we focus on Case (D-2), where two Qs = I reso- 
nant surfaces are located a small distance apart. The 
results are similar in other cases with different Qs, pro- 
vided that TOpeak IS Several times larger than the lowest 
possible mode number m. The evolution of the individ- 
ual Fourier modes is described in terms of the kinetic and 
magnetic energies of their perturbation components, 



„j 100 200 300 400 500 600 700 800 900 1000 



£^r„ = |V0™,„|^ and EZ^, = mm,n\' 
and the corresponding nonlinear growth rates. 



(15) 



kin ^ dl^-^m"n 
rrn.nK ) 2dt 



dlnK^J^s 
and ^llfAt) = 2^;"-" . (16) 



(these are amplitude growth rates, hence the factor 1/2). 
In Eq. (HSl), |/„,,„|2 = /g^dr r C„|/„,,„(r)|2, with 
Cra=Q = 47r and C,„^o = Stt. 

In Fig. [5] the evolution of (a) the energies and (b) 
growth rates of the m = 1 and m = modes in Case 
(D-2) is shown. Note that the m = mode considered 
here measures only the profile perturbation, excluding the 
equilibrium profile. During phase (i) the linear mode 
structure of the m = 1 mode is gradually established. In 
the present case, this process is not fully completed by 
the time the nonlinear drive (ii) sets in, as can be inferred 
from the fact that and 7™i^ are not equal and still 
vary in time. 

The nonlinear driving phase (ii) begins when the 
fastest growing modes (not shown) reach sufficiently large 
amplitudes so that they start to drive slower modes 
through nonlinear coupling. For instance, a typical driv- 
ing term in the convectivc nonlinearity 0] in Eq. ([T|) 
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FIG. 8: (Color online). Evolution of the m = 1 and m = 
modes in Case (D-2) . In (a) the evolution of the magnetic and 
kinetic energies E'^''<^ and S'''" [Eq. JlS])] is shown. In (b) the 
magnetic and kinetic growth rates 7™"^ and 7*^'" [Eq. (|16|l ] 
are plotted. The main stages relevant to the evolution of the 
771 = 1 mode are: (i) establishing the linear mode structure 
and linear growth, (ii) nonlinearly driven growth, (iii) recon- 
nection in the inter-resonance region (annular collapse). In 
diagram (b), the dashed horizontal line in phase (i) indicates 
the linear growth rate 7iin('Ti = 1) = 4.0 x 10~^ from Fig. [2] 
The dashed horizontal line during stage (ii) indicates the ex- 
pected growth rate due to nonlinear driving. Sup = 10*^, 
Re-Yip = 10*, ^'o,m>o = 10~^, and the initial phases for modes 
m = 6-10 are i9om = {0, 0, tt, tt, 0}. Similar behavior is also 
found for cases with qs > 1, such as (D-1). 



is 

^'^{r) sin(mat?*)e''"* x Mr) sm{mb'd^)e'^''* (17) 
= (?/'>6/2) [cos (m-?9,) - cos (m+?9,)] e('^"+''')*, 

where rria and m?, are the two driving modes with linear 
growth rates 7a and jt, and = ma^m}, are the mode 
numbers of driven modes. The growth rate of the driving 
term is 7drive(m='=) = 7^ -f 7fc- 

In Fig. [8] the m = 1 mode switches between (al- 
most) linear and nonlinearly driven exponential growth 
around t w 300. Here, the driving is primarily due to 
the coupling between the modes m = 8 and m = 7. 
Consequently, according to Eq. ([TT]), the expected en- 
hanced growth rate is 7drivc("i = 1) « 7drivc("i = 0) = 
27iin("^ = 8) = 17.5 X lO^-^. The level of nonlinear driv- 
ing can conveniently be inferred from the growth rate of 
the m — mode (only magnetic energy). The m — 
mode is not an unstable eigenmodc and its evolution is 
entirely a result of nonlinear driving by higher-m modes. 

The annular collapse, labeled by (iii) in Fig. [51 begins 
already during the driving phase and continues beyond 
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FIG. 9: (Color online). Evolution of the growth rates of the 
m = 1 and m = modes in Case (D-2). As in Fig. [8jb), but 
with different initial conditions: (a) The perturbation ampli- 
tude is V&o = 10^^^ for all modes, and the initial phases for 
the five dominant modes m = 6-10 are i?om = {0, 0, tt, tt, 0}. 
(b) The perturbation amplitude is "fo = lO"^'^ for all modes, 
and the initial phases for the dominant modes m — 6-10 are 
i^Om ~ {0,0,0,0,0}. These results were obtained with a re- 
duced number of modes, m — 0-31, and are not valid far 
beyond t ~ 1700 rup. Sup = 10^ Reup = 10**. 



it. The dynamics during this stage were described above 
in Section IV^ 

In Fig. [8] neither the m — nor the m ~ 1 mode reach 
the expected growth rate, 7drivc = 17.5 x lO^'^. The 
reason for this Hes in the initial perturbation. We have 
determined two factors that need to be considered: (A) 
the perturbation amplitude, and (B) the phase relations 
between the driving modes. 

The effect of (A) can be seen by comparing Fig. [Hl^b) 
with Fig. IHlJa). Figure [9] shows results for Case (D-2) 
that were obtained with a lower perturbation amplitude 
= 10-12 [ij^ Fig. H = 10"'^]. Now, there is more 
time for all modes to establish the linear mode structures. 
Clearly, the nonlinear driving in Fig.[21[a) reaches a higher 
level than in Fig. [Hb). 

Factor (B) implies that the growth rate 7drivc may also 
stay below 27pcak when there are several higher-m modes 
with growth rates approximately equal to 7pcak and when 
the phase relations between these modes are "unfavor- 
able." Thus, 7drive Can be increased by "aligning" the 
phases of the fastest growing modes in the spectrum. The 
result of aligning the five modes m = 6-10 is shown in 
Fig. [DJb) . Clearly, the driving of the to = 1 mode is now 
much more effective. 

Strictly speaking, the effect of aligning the relative 
phases is merely to reduce the time needed to establish 
the nonlinear driving. Its effect on the growth rate is 
only temporary. Note that the growth of the to = 




mode structure durincj 
nonlinear driving ; 



FIG. 10: (Color online). Effect of the nonlinear driving on 
the m — 1 mode structure in Case (P-2) [cf. Fig. IHb)]. The 
two fastest growing modes are m = 8 and m = 7 [cf. Fig. [2]. 

In (a) the linear eigenmode structure V'l^i is shown. In (b) 
and (c) the typical structure of the driving terms can be seen. 
Here, the terms {drit>H,8)[T / r)4>7 ,7 and (8/r)?/)8,8(9r(^7,7) are 
shown, representatively for the convective nonlinearity [iji^fji] 
in Eq. ([T}. In (d) the mode structure during the nonlinear 
driving phase is plotted. 



mode is not affected by phase relations between the driv- 
ing modes, as is to be expected. 

The effect of the nonlinear driving on the to = 1 mode 
structures can be observed in Fig. 1101 Note in partic- 
ular that despite the radial localization of the driving 
terms [Fig. [TU^b) and (c)], the to = 1 mode [Fig. [TOlJa)] 
as a whole grows at an enhanced rate. Once the driving 
is established the mode structure varies only minutely 
[Fig. llOf d)]. Thus, a new nonlinear mode structure is 
formed. One consequence of this global effect of the lo- 
calized driving is that in cases where the driving term 
happens to have the opposite sign, it induces a switching 
of the sign of the global m ~ 1 mode structure at the 
onset of the nonlinear driving phase. Such an event can 
be observed in Fig.[9jb), where S"™!^ performs an under- 
and overshoot as a result of sign reversal. 



VI. DISCUSSION AND CONCLUSIONS 

In this paper we have studied the linear instability and 
nonlinear dynamics of resistive DTMs for small inter- 
resonance distance D12 . The results may be summarized 
as follows. 

The linear growth rates 7iin(n^) of 1ow-to modes were 
found to increase with to. The power law 7iin oc m?/^ 
predicted by the linear theory developed in Ref. [l3| 
agrees with numerical results. The linear growth rates 
of high-TO modes decrease with to. These modes seem 
to be outside the scope of the existing theory. A semi- 
empirical formula for the mode number of the fastest 
growing mode TOpoak was proposed [Eq. (jl3|) ]. based on 
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results of Ref. [3] ■ The estimate for mpeak is valid only 
for resistive DTMs. Nevertheless, broad spectra and high 
values mpoak are also found when reconnection is medi- 
ated by electron inertia instead of resistivity, as is shown 
in a separate paper [2lj. Also, Eq. (fT3|) does not ap- 
ply to cases with more than two resonant surfaces. For 
instance, for triple tearing modes (TTMs) both mpoak 
and the number of unstable modes arc larger than for 
DTMs under comparable conditions, indicating that a 
TTM tends to be an even stronger instability [23 • 

The dominance of several high-m DTMs in config- 
urations with small D12 leads to an annular collapse 
with small island structures. This turbulent annular col- 
lapse shows that the nonlinear DTM dynamics in cases 
with small inter-resonance distances (dominant high-m 
modes) are different from those studied previously by 
other authors using large inter-resonancc distances (dom- 
inant lowest-m mode) . Our results may be of interest for 
scenarios, where (jmin is gradually lowered by increasing 
the off-axis current. Experimentally, an apparently qui- 
escent passage through low-order resonant surfaces and 
thus access to regimes with larger inter-resonance dis- 
tance, was found to be possible by applying additional 
external drive (e.g., Ref. Q). MHD activity was re- 
ported after the inter-resonance distance had grown sub- 
stantially and the dynamics seem to be dominated by 
low-m modes. It is possible that earlier high-m activity 
predicted by our simulations is cither prevented by the 
external drive (not included in our model) or has escaped 
detection. Through recent progress in plasma diagnos- 
tics, which allows to detect high-m magnetic islands on 
low-order resonant surfaces [23|, experimental checks of 
our simulation results may be feasible in the near future. 

Before entering the fully nonlinear regime, the fastest 
growing modes drive slower modes. For Qs = 1, this in- 
cludes the global m = 1 resistive internal kink mode, de- 
spite the fact that the driving is radially localized between 
the resonant surfaces. This nonlinear driving experienced 
by slower modes is of particular interest for the m = 1 
internal kink mode. The nonlinear driving implies that 
the m = 1 mode appears significantly earlier than would 
be expected from its linear growth rate. Hence, the fast 
growing high-m DTMs effectively provide a trigger mech- 
anism for the m = 1 mode, similarly to the TTM case 
described in detail in Ref. [la|. Note that, by the time 
the m = 1 mode becomes observable, the inter-resonance 
region has already undergone reconnection. The lack of 
magnetic shear allows the m = 1 mode to grow as an 
ideal internal kink instability. In contrast to the TTM 



case, there is no sawtooth crash associated with the ac- 
tivity of the m = 1 mode when the inter-resonance dis- 
tance is small. However, the rapid excitation of m = 1 
oscillations observed in our simulations may be useful in 
experiments to determine the instant in time when Qnmi 
passes through the qs = 1 resonant surface. 

It is important to note that the instability of a broad 
spectrum of modes implies that the details of the non- 
linear dynamics depend on the initial conditions used in 
the simulation. For instance, the efficiency of the non- 
linear driving was shown to depend on the phase re- 
lations between the fastest growing modes (cf. Fig. [9]). 
Furthermore, the phase relations determine where the 
first magnetic islands form. In principle, it is possible 
to produce the first magnetic islands in a poloidally lo- 
calized region. Although, this discussion of initial per- 
turbations may seem rather academic, we believe that it 
may bear practical importance. For instance, pellet or 
neutral beam injection may provide a localized magnetic 
perturbation and possibly form robust helical structures 
such as "snakes" [2J]. Furthermore, if the DTMs are ex- 
cited by micro-turbulence, modes with m > mpoak may 
dominate the annular collapse, in contrast to the present 
study, where all modes were perturbed with similar am- 
plitude so that mpoak became dominant. 

In conclusion, the results presented in this paper may 
be relevant to the understanding of the MHD activity 
near qmin and may thus be of interest for studies on stabil- 
ity and confinement of advanced tokamaks. The results 
motivate further research in this direction with more re- 
alistic models. 
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